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Abstract — This paper presents three feedback controllers 
that achieve an asymptotically stable, periodic, and fast 
walking gait for a 3D bipedal robot consisting of a torso, 
revolute knees, and passive (unactuated) point feet. The 
walking surface is assumed to be rigid and flat; the contact 
between the robot and the walking surface is assumed to 
inhibit yaw rotation. The studied robot has 8 DOF in the 
single support phase and 6 actuators. In addition to the re- 
duced number of actuators, the interest of studying robots 
with point feet is that the feedback control solution must 
explicitly account for the robot's natural dynamics in order 
to achieve balance while walking. We use an extension of 
the method of virtual constraints and hybrid zero dynamics, 
a very successful method for planar bipeds, in order to si- 
multaneously compute a periodic orbit and an autonomous 
feedback controller that realizes the orbit, for a 3D (spa- 
tial) bipedal walking robot. This method allows the com- 
putations for the controller design and the periodic orbit to 
be carried out on a 2-DOF subsystem of the 8-DOF robot 
model. The stability of the walking gait under closed-loop 
control is evaluated with the linearization of the restricted 
Poincare map of the hybrid zero dynamics. Most periodic 
walking gaits for this robot are unstable when the controlled 
outputs are selected to be the actuated coordinates. Three 
strategies are explored to produce stable walking. The first 
strategy consists of imposing a stability condition during the 
search of a periodic gait by optimization. The second strat- 
egy uses an event-based controller to modify the eigenvalues 
of the (linearized) Poincare map. In the third approach, the 
effect of output selection on the zero dynamics is discussed 
and a pertinent choice of outputs is proposed, leading to 
stabilization without the use of a supplemental event-based 
controller. 



I. Introduction 

The primary objective of this papeiQ is to contribute to 
the feedback control of 3D bipedal robots that do not rely 
on large feet and slow movement for achieving stability of a 
walking gait. We assume here an unactuated point contact 
at the leg end and, for a simple 5-link robot, seek a time- 
invariant feedback controller that creates an exponentially 
stable, periodic walking motion. Our approach is based on 
an extension of the method of virtual constraints, which 
was developed in [2], [3], [4], [5] for planar robots, and is 
extended here to the case of spatial robots. Virtual con- 
straints are holonomic constraints on the robot's configu- 
ration that are asymptotically achieved through the action 
of a feedback controller. Their function is to coordinate 
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the evolution of the various links of the robot throughout a 
stride — which is another way of saying that they reduce the 
degrees of freedom. By using virtual constraints to achieve 
link coordination on a bipedal robot, different gaits can be 
more easily programmed than if the links were coordinated 
by hardware constraints. 

The work most closely related to ours is [B], where the 
control of a 3D walker was decomposed into the study 
of its motion in the sagittal plane and the frontal plane; 
see also [7] for a related decomposition result on control 
in the frontal plane. The method of virtual constraints 
was applied in [5J to regulate the sagittal plane motion 
of the biped, while an inverted pendulum approximation 
of the dynamics was used to design a controller for the 
frontal plane. An event-based controller was then intro- 
duced to synchronize the phasing of the independently de- 
signed sagittal and frontal plane controllers. The overall 
closed-loop system was shown to be stable through simu- 
lation and subsequently through experimentation. In our 
approach, we do not decompose the model into sagittal 
and frontal plane motions, and coupling of the sagittal 
and frontal plane dynamics is introduced into the controller 
from the very beginning. 

A very interesting study of the feedback control of under- 
actuated spatial robots has been given in [5] , where a con- 
troller for a five-link 3D robot with unactuated point feet 
has been designed on the basis of linearizing the robot's dy- 
namic model along a periodic orbit. So that the controller 
would be time- invariant, the orbit was parameterized with 
a configuration variable that is strictly monotonic through- 
out a normal gait, as in [2], [3], [4], [5], before linearization 
was applied. The (within-stride) control law is designed on 
the basis of a discrete-time approximation of the linearized 
model, which makes stability of the closed-loop system dif- 
ficult to assess. 

Other important work includes [5] and references therein, 
where the analysis of passive spatial bipeds is presented. 
The emphasis in their work is on energy efficiency and un- 
deractuation; the role of feedback control in achieving a 
wide range of behaviors is not emphasized. On the other 
hand, the work in [TU], [TT] seeks energy efficiency and a 
large basin of attraction under the assumption of full ac- 
tuation; in particular, full actuation between the leg and 
ground is assumed (pitch, roll and yaw) , as opposed to the 
unactuated assumption made here. Very careful stability 
analysis of the closed-loop system is provided through ge- 
ometric (Routhian) reduction. This work is taken one step 
further in [12], where, starting from a 2D (sagittal plane) 
passive limit cycle, the authors use geometric reduction to 
first achieve control of the frontal plane motion and then 
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a second stage of geometric reduction to achieve steering 
within the walking surface. 

To the best of our knowledge, other work on the con- 
trol of spatial robots either assumes full actuation or does 
not provide significant analysis of the closed-loop system. 
There are many control strategies based on the zero mo- 
ment point ZMP [13] . with one of the more famous users 
being the robot ASIMO [14]. In this approach, a desired 
trajectory of the ZMP is defined and successive inner con- 
trol loops are closed on the basis of the ZMP. In the work 
of |15j , predictive control is performed on the basis of the 
position of the center of mass and a simplified model of the 
robot in order to achieve a desired ZMP trajectory. Re- 
cently, on-line adjustment of the ZMP has been added [15] : 
this control method is implemented on the robot HRP2. 
The control of the ZMP ensures that the supporting foot 
will not rotate about its extremities, but this does not en- 
sure stability in the sense of convergence toward a periodic 
motion, as proved in |17) . 

II. Model 

A simplified model of a spatial bipedal robot is given 
here. The model was chosen to be complex enough to cap- 
ture interesting features of gait control that do not occur 
in planar robots, and simple enough that the presentation 
of the ideas will remain transparent. It is our expecta- 
tion that the ideas presented here apply to a wider class of 
bipeds, but proving such a conjecture is not the objective 
of this paper. 

A . Description of the robot and the walking gait 

The 3D bipedal robot discussed in this work is depicted 
in Fig- [D It consists of five links: a torso and two legs with 
revolute one DOF knees that are independently actuated 
and terminated with "point-feet" . Each hip consists of a 
revolute joint with two degrees of freedom and each degree 
of freedom (DOF) is independently actuated. The width 
of the hips is nonzero. The stance leg is assumed to act as 
a passive pivot in the sagittal and frontal planes, with no 
rotation about the z-axis (i.e., no yaw motion), so the leg 
end is modeled as a point contact with two DOF and no 
actuation. As discussed below, this model corresponds to 
the limiting case of robot with feet when the size of the feet 
decreases to zero. The unactuated DOF at the leg ends 
correspond to the classical DOF of an ankle. The DOF 
corresponding to the swing-leg ankle are not modelled. In 
total, the biped in the single support phase has eight DOF, 
and there are two degrees of underactuation. 

In a more complete model with feet, one would include 
at a minimum two degrees of freedom in the ankles, cor- 
responding to motion in the sagittal and frontal planes. 
Moreover, it would be assumed that the friction between 
the foot and the ground is sufficient to prevent sliding, and 
hence in particular, rotation of the foot in the yaw direction 
[18] . [T5] . In most control studies, which assume flat-footed 
walking, the key limiting factor is the determination of an- 
kle torques that respect the ZMP condition, namely the 
ground reaction forces must remain with the convex hull 



of the foot [20], [2T], [13] . The lighter and smaller feet, the 
tighter are the constraints on the allowable ankle torques, 
increasing the difficulty of determining feasible walking tra- 
jectories and subsequently, stabilizing feedback controllers. 

The objective of our study is two-fold: to show that sta- 
ble flat-footed walking is possible with zero ankle torques, 
thereby removing an important obstacle to previous stud- 
ies on walking; and to prepare for a future study of walking 
gaits that allow phases with rotation about the heel and toe 
[22], [17]. In the limiting case of a point foot, the allowable 
ankle torque becomes zero, leading to an unactuated con- 
tact at the leg end. Given that the envisioned application 
of our results is to robots with feet, where yaw rotation is 
naturally inhibited by friction, it makes sense to assume in 
the point-foot model that yaw rotation is not allowed. 

In summary, the following assumptions are made in the 
present study: 

« Each link is rigid and has mass. 

• Walking consists of two alternating phases of motion: 
single support and double support. 

• The double support phase is instantaneous and occurs 
when the swing leg impacts the ground. 

• At impact, the swing leg neither slips nor rebounds. 

• The swing and stance legs exchange their roles at each 
impact. 

• The gait is symmetric in steady state. 

• Walking takes place on a flat surface. 

A more detailed list of hypotheses is given in [23J Chap. 3] 
for planar robots, and with the obvious modifications for 
spatial robots, those hypotheses apply equally well here. 




q 4 



Fig. 1 

a five-link 3d biped with point feet in support on leg-1. 
There is no yaw motion about the stance leg end and the 
DOF at the leg end are unactuated. For simplicity, each 

LINK IS MODELLED BY A POINT MASS AT ITS CENTER. 

Since the gait is composed primarily of single support 
phases, the variables used to describe the robot are adapted 
to this phase of motion. The robot is represented as a tree 
structure. The stance foot, which is fixed on the ground, is 
the base of the tree structure. A set of generalized coordi- 
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nates q = [qi, . . . , q&]' is shown in Fig. [TJ Absolute angles 
((71,(72) are roll and pitch angles of the stance leg, respec- 
tively. Angles (73 and q% are the relative joint angles of the 
stance-leg knee and swing-leg knee, respectively. Angles (74 
and (75 are the joint angles of the stance leg relative to the 
torso along the y-axis and the x-axis, respectively, and an- 
gles qe and (77 are the joint angles of the swing leg relative 
to the torso along the x-axis and the y-axis, respectively. 
The coordinates (q\, (72) are unactuated (due to the passive 
contact), while (53, . . . ,qs) are independently actuated. 

The position of the robot with respect to an iner- 
tial frame is defined by adding the four variables q e = 
[q / ,x st ,y s t,z st ,qo, s t\', where x st , y s *and z st are the Carte- 
sian coordinates of the stance foolQ and (70. s t defines the 
rotation along the z-axis of the stance leg. These variables 
are constant during each single support phase. 

We have chosen to define the generalized coordinates 
with respect to the contact point of the current stance foot. 
When leg-2 is the supporting leg, the variables are defined 
as shown in Fig. [5] and the same notation is employed as 
when the supporting leg is leg-1, viz. Fig. [TJ Hence, at 
each leg exchange (i.e., impact), the variables q e undergo a 
jump due to the change in location of the reference frame. 

During each single support phase, only one set of coor- 
dinates is used, depending on which leg is the supporting- 
leg. In double support, either set of coordinates may be 
used. The transformation from one set of coordinates to 
the other is nonlinear [8], but it can be computed in closed 
form by standard meanqj. 




Fig. 2 

A FIVE-LINK 3D POINT-FEET BIPED IN SUPPORT ON LEG-2. THE DOF 
AT THE LEG END (FOOT) ARE NOT ACTUATED. 



The legs exchange roles from one step to the next. If T 
is the duration of a step, on a periodic walking cycle, due 

2 The leg ends are referred to as feet or point feet. 

3 The transformation from one set of coordinates at the end of a 
step, for example, to the other set of coordinates is done as follows. 
Compute the orientation and the angular velocity of the swing leg 
shin. From this, one deduces qg, qi and 92 that are compatible with 
this orientation, and then one deduces 90, 91 and 92 yielding the 
angular velocity of the swing shin. The angles 93 to qg exchange 
their roles viz [93,94,95,56,97,98] -> [98, 97, 96, 95, 94, 93] ■ 



-30 - 




X 

Fig. 3 

Trajectory of the robot's center of mass projected onto 
the 2i-2/-plane, showing that a net rotation is possible even in 

a model without explicit yaw rotation. this motion arises 
from rotational coupling between the frontal and sagittal 

PLANES. 



to the choice of coordinates in Figs. [JJand[21 we must have 



qi(t + T) = -q x {t) 
q 3 (t + T) = q 3 (t) 
q 5 (t + T) = -q 5 (t) 
q 7 (t + T) = q 7 (t) 
and 

qo,st(t + T) = -q 0tSt (t). 



q 2 (t + T) = q 2 (t) 
q 4 {t + T) = q 4 (t) 
q 6 {t + T) = -q 6 (t) 
q s (t + T) = q 8 (t) 



(1) 



The last condition yields a motion along the x-axis. 
Remark: Even though the model does not include rota- 
tion about the supporting foot, nor any other vertical axis, 
such as at the hip, coupling between the rotations in the 
sagittal and frontal planes can yield a net rotation about 
the vertical axis from one step to the next. Hence it is 
important to keep track of (70, s *. As an example, Fig. [3] 
shows the projection onto the x-y-plane of a solution of 
the model yielding a circular motion. This was obtained 
by modifying the control law of I Villi to have unequal step 
lengths on the right and left legs. 

B. Dynamic model 

The dynamic models for single support and impact (i.e., 
double support) are derived here assuming support on leg 
1. The models for support on leg 2 can be written in a 
similar way by using a hip width of — W in place of W 
(this places the swing leg on correct side of the stance leg). 
The Euler-Lagrange equations yield the dynamic model for 
the robot in the single support phase as 



D(q)q + H(q,q)=Bu 



O2X6 
-^6x6 



U, 



(2) 



where D(q) is the positive-definite (8 x 8) mass-inertia ma- 
trix, H(q,q) is the (8 x 1) vector of Coriolis and gravity 
terms, B is an (8 x 6) full-rank, constant matrix indicat- 
ing whether a joint is actuated or not, and u is the (6 x 1) 
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vector of input torques. Following standard practice in the 
literature, the double support phase is assumed to be in- 
stantaneous. However, it actually consists of two distinct 
subphases: the impact, during which a rigid impact takes 
place between the swing foot and the ground, and coordi- 
nate relabeling. During the impact, the biped's configu- 
ration variables do not change, but the generalized veloci- 
ties undergo a jump. The derivation of the impact model 
in double support phase requires the use of the vector q e . 
Conservation of angular momentum of the robot about the 
end of the swing leg during the impact process, in combi- 
nation with the swing leg neither slipping nor rebounding 
at impact, yields 



(3) 



where q~ and q+ are the extended velocities before and 
after the impact, respectively, F sw is the reaction force at 
the contact point, D e is the extend mass-inertia matrix, 
and E sw = [x sw ,y sw ,z sw ,q 0tSW ]' is the Jacobian matrix 
for the position of the swing foot and its orientation in the 
x — y-plane. Analogously to [2], the overall impact model 
is written as 

q + = A,(g") (4) 

and 

q+ = A 4 (q-,q-), (5) 

and is obtained from solving <j3j> and projecting down to 
the generalized coordinates for support on leg 2. 



Zsw{q) = & Xaw(q) > 



\4t " 




D e -E' sw 


-1 


' D e q- ' 


± sw 




E S w O4 x 4 




04xi 



Define state variables as 



and let 



(I 



and x- 



, where the subscript j € {1,2} 



denotes the stance leg number. Then a complete walking 
motion of the robot can be expressed as a nonlinear system 
with impulse effects, as shown in Fig. |4] and written as 



fi(xi) + gi{xi)ui Xi <£ Si 
-- Ai(a;jf ) x\ G Si 

h{%2) + #2(2:2) U 2 X2 £ S 2 



(6) 



x+ = A 2 (x 2 ) 



G S 2 



where S± = {(q,q)\z sw (q) = 0, x 8W (q) > 0} is the switch- 
ing surface, 



AO*) 

and 



q 



-D- 1 (q)H(q,q) 



x+ = A 1 (x 1 ) 



A,(g-) 

A 4 (q-,q-) 





D-\q)B 



When leg-2 is the support leg, the same derivation produces 
S%, h, 92 and A 2 . 

III. Virtual constraints 

The method of virtual constraints, which has proven very 
successful in designing feedback controllers for stable walk- 
ing in planar bipeds [2], [3], [4], [5], will be applied to the 




x 2 = A l( x l ) 
"1 ( ±2 = h ( x 2) + g2(x 2 )u 2 ) 

x+ = A 2 (x~) 



Zsw(q) = & x sw (q) > 
Fig. 4 

Bipedal robot's dynamic model as a hybrid system. 



3D biped of the previous section. In this method, one holo- 
nomic constraint per actuator is proposed in the form of an 
output that, when zeroed by a feedback controller, enforces 
the constraint. The most direct form of the constraint is 



V = HQ) =Qa- h d (0), 



(7) 



where q a = [q3,q4,,q5,qe,Q7,Qs]' is the vector of actuated 
coordinates, 9 = 0(q) is a quantity that is strictly mono- 
tonic (i.e., strictly increasing or decreasing) along a typical 
walking gait, and hd{9) is the desired evolution of the ac- 
tuated variables as a function of 9. Roughly speaking, 9 is 
used to replace time in parameterizing a periodic motion of 
the biped. In a forward walking motion, the x-coordinate 
of the hip is monotonically increasing. Hence, if the virtual 
stance leg is defined by the line that connects the stance 
foot to the stance hip, the the angle of this leg in the sagit- 
tal plane is monotonic. When the shin and the thigh have 
the same length, the angle of the virtual leg in the sagittal 
plane can be selected as 



-<?2 - <73/2 



(8) 



(the minus sign is used to make 9 strictly increasing over a 
step). 

The torque u* required to remain on the virtual con- 
straint surface corresponding to q a — hd{9) can be com- 
puted afl 



d9 2 



dq 



(9) 

This leads to an input-output linearizing controller to 
asymptotically drive the state of the robot to the constraint 
surface, assuming it does not initially start there [5S] [531 
Chap. 5], 



u = u 



which results in 



dh 



^d-'b) C-^y + —y), 



.. - ~ ., . , , (10) 

oq ) e £ 



.. K d K p 

y H v + —ry = 0. 



(11) 



In other words, determining the constraints is equivalent 
to the design of a feedback controller in the single support 

4 As shown in [24] . [231 pp. 60], an expression without inversion of 
the (8 X 8) mass-inertia matrix D is also possible. The expression 
given in (01 is more compact. 
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phase, up to the choice of the gains K p > 0, Kd > 0, and 
e > such that (fTTjl is exponentially stable and converges 
sufficiently rapidly with respect to the duration of a single 
support phase; see [23l Chap. 4]. 

The next objective is to determine the behavior of the 
robot under the virtual constraints. This task is simplified 
by noting that enforcing the virtual constraints, y = h(q) = 
0, results in q a = hd(9) and reduces the dimension of the 
dynamics. 

Let q u = [qi,9]' denote the unactuated joints and q a = 
[937 94j 95: 1e, 17, <7s]' denote the controlled joints, which are 
selected here to be the actuated joints. A linear relation 
exist between q, q u and q a , 



q = T 



Qv 



(12) 



where T is an (8 x 8) invertible matrix. Then ([2]) can be 
rewritten as 



T'D(q)T 



Qa 



-T'H(q, q) = T'Bu = 



02x6 
^6x6 



u, (13) 



The first two lines of the RHS of this equation are zero, 
yielding 



Dn{q)q u + D 12 {q)q a + Hi(q,q) = 2x i, 



(14) 



where D\\ is the (2 x 2) upper left sub-matrix of T'D(q)T, 
D12 is the (2 x 6) upper right sub-matrix of T'D(q)T and 
Hi(q, q) consists of the first two lines of T'H(q, q). Substi- 
tuting the expressions of q a , q a and q a corresponding to the 
virtual constraints, the dynamic model of the single sup- 
port phase is now reduced to a low-dimensional, 2-DOF, 
autonomous system, 



Du{q u ) 



+ D 12 (qu)(^t 
Hi(q u ,q u ) = 0, 



d 2 h d hi 



(15) 



which is called the swing phase zero dynamics [26| 
Chap. 5]. 

One can clearly see that the dynamic properties of the 
swing phase zero dynamics depend on the particular choice 
of the virtual constraint y = q a — hd(0) = 0. How to deter- 
mine a choice for hd(9) that results in a periodic walking 
motion is summarized in the next section. 

IV. Design for a Symmetric Periodic Gait 

The objective of this section is to design virtual con- 
straints q a = hd(0) that correspond to a periodic motion of 
the robot. The gait considered is composed of single sup- 
port phases separated by impacts as described in Fig. 
The legs exchange roles from one step to next, and due to 
symmetry, the study of a gait can be limited to a single 
step and the use of the symmetry relation ([1]). 

A . Virtual constraints and Bezier polynomials 

The problem of designing the virtual constraints will be 
transformed into a parameter optimization problem as in 



[23l Chap. 6]. Here, our main goal is to obtain a peri- 
odic motion; optimality is not so crucial. To simplify the 
optimization process, the number of variables used in the 
optimization problem is first reduced. This is accomplished 
by exploiting boundary conditions that arise from periodic- 
ity. Bezier polynomials are parametric functions that allow 
one to easily take into account boundary conditions on the 
configuration and velocity at the beginning and end of a 
step. 

The initial and final configuration and velocity of the 
robot for a single support phase are important for defining 
the passage between the single and double support phases. 
Because the terminal configuration of the robot is chosen 
to be the instant before the double support configuration, 
both legs are in contact with the ground and therefore 
only seven independent variables are needed to describe 
this configuration (a closed kinematic chain). These vari- 
ables parameterize the final configuration of the first step 
denoted qf. The eight joint velocities qf are independent 
and are also added. 

Knowing the final state of the single support phase, the 
impact model Q and (0 determines the initial state of 
the ensuing single support phase. The symmetry condition 
(HJ) then gives the initial state of the first step: qi, q\. The 
initial orientation (qi)o,st of the robot is calculated such 
that the orientation for the second step is symmetric to the 
orientation for the first step in order that no yaw rotation 
is observed during the nominal (periodic) gait. 

To obtain a periodic gait, the single support must be 
such that the state of the robot evolves from qi, qi to qf, 
qf. For given desired initial and final state values, virtual 
constraints can be easily deduced to connect the desired 
values of the actuated variables. However, the evolution 
of the unactuated variables is known only by integration of 
the dynamics (|15[) ; a desirable dynamic behavior is imposed 
on these variables by the use of equality and inequality 
constraints in the optimization process. 



23| B. Specifics 



Here, Bezier polynomials of degree 3 are chosen to de- 
fine the virtual constraint^]. The virtual constraints are 
expressed as functions of the variable 9; see ((SJ. From qi 
and qf, the initial and final values of 9, denoted 9i and 9 y, 
can be calculated. Let 



3! 



fe!(3 - fc)! 



3-k 



where 



(16) 



(17) 



is the normalized independent variable. The coefficients 
of the Bezier polynomials, ctk, are (6 x 1) vectors of real 
numbers. They must be determined so as to join (qi) a 
to (qf) a and (qi) a to (qf) a , (the additional subscript "a" 

5 A degree greater than 3 can also be chosen, in which case the 
number of optimization variables increases I27| . 
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denotes the actuated variables) when 6 varies from 9i to 
Of, yielding 

a = h d (6i) = (qi) a 

«i = (Qi)a + = (*)- + 6 -^ {Jt t , 1B , 

a* = (<lf)a-^^(6f) = ( q f) a -^^ (18) 

«3 = h d {6f) = (q f ) a - 

The evolution of the unactuated variables is calculated by 
integration of the dynamic subsystem (fl"5j) . that is, the 
stance phase zero dynamics, starting from the initial state 
{li)u = [{li)i-,0i\' and terminating at = Of, where (qi)i 
denotes the initial value of q\ for the first step. 

When the evolution of the unactuated variables is cal- 
culated, because the evolution of the actuated variable is 
given by (0 and (|16j) . the torque required to zero the con- 
straints, i.e, u* in ©, can be calculated by lines 3 through 
8 of (fl"3|) . and the ground reaction force F st expressed in 
the inertial reference frame (see Fig. [IJ can be calculated 
as well. 

The search for a periodic walking motion can now be 
cast as a constrained nonlinear optimization problem: Find 
the 15 optimization parameters prescribing (qf,qf) that 
minimize the integral-squared torque per step length^, 

1 f T > 

J = — u* u*dt, (19) 
L Jo 

where T is the walking period and L is the step length, 
while satisfying symmetry ([1]), and subject to the following: 
inequality constraints 

• 6 is strictly increasing (i.e, > along the solution) ; 

• the swing foot is positioned above the ground (z sw > 0); 

• a no-take-off constraint, F st (3) > 0; 

• a friction constraint, y/F st (l) 2 + F st (2) 2 < /j,F st (3); 
equality constraints 

and a set of conditions imposing periodicity, 

H ( T ) = (<?/)i 
4i(T) = (<?/)! 
0\T) = Of, 

where qi(t) and 0(t) result from the integration of the zero 
dynamics and the walking period T is such that 0(T) — Of- 
The above procedure can be performed in MATLAB 
with the FMINCON function of the optimization toolbox. A 
fixed-point solution x* — [q*^q*j]' minimizing J defines a 
desired periodic walking cycle (or nominal orbit). The cri- 
terion being optimized (|19p has many local minima and the 
optimization technique used is local. Thus, the obtained 
optimal periodic motion depends on the initial set of opti- 
mization parameters. 

6 Torquc being proportional to current in a DC motor, integral- 
squared torque is a rough approximation of energy dissipated in the 
motors. 



C. An example periodic motion minimizing integral- 
squared torque 

The physical parameters of the 3D biped studied here are 
given in TablelU For these parameters, a periodic orbit was 
computed following the technique presented in the previous 
subsection. We obtained a periodic motion defined by x* = 
(<lf,Qf), where 

q* f = [-0.0174,-0.34038,0.3820,-0.2940,0.0602,0.0487, 

-0.5077,0.1688]', 
q* f = [-0.4759,-1.1825,0.0997,0.2785,-0.1000,0.1000, 
1.398,0]'. 



g 


W 


LI 


L2 


L3 


ml 


m2 


m3 


9.81 


0.15 


0.275 


0.275 


0.05 


0.875 


0.875 


5.5 



TABLE I 

Parameters for the 3D bipedal robot (in MKS). 




Fig. 5 

Stick-diagram of the optimal motion with respect to the 
torque criterion. 

A stick-figure diagram for the first step of the periodic 
walking gait is presented in Fig. [S] The walking gait has a 
period of T = 0.39 seconds, a step size of L = 0.176m, and 
an average walking speed of 0.447 m/sec, or 0.745 body 
lengths per second. The step width is 0.156m, close to the 
hip width. The nominal gait's joint profiles and angular 
velocities over two consecutive steps are shown in Figs. [5] 
and[7l respectively. The unactuated and actuated variables 
are presented; note that 9 is monotonic over each step. 
Fig. |8] shows the torque required to produce the periodic 
motion, which is less than lONm for each joint. Fig. [9] 
shows the profile of the ground reaction force on the stance 
foot and the profile of the swing leg tip; this figure shows 
that the inequality constraints are satisfied on the nominal 
motion. 

V. Evaluating the Stability of a Walking Cycle 

The stability of a fixed-point x* can be tested numeri- 
cally by linearizing the Poincare map about the fixed-point 
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time (s) 



time (s) 



Fig. 7 



Joint rate profiles of the obtained periodic motion over 
two steps, where the small circles represent qj. 




time (s) 



time (s) 



Fig. 8 

Torques profiles of the obtained periodic motion over two 

STEPS 



as presented in [I] . This numerical stability test has a high 
computational cost, however, because it requires the esti- 
mation of the Jacobian of the Poincare map, in a space of 
dimension 2n-l, where n is the number of independent joint 
coordinates; here n = 8. We propose a slight modification 
of the control law in order to be able to study the stability 
of the closed-loop system in a reduced-dimensional state 
space. 

A. Hybrid zero dynamics (HZD) and a stability test in a 
reduced space 

The control law (| 10[) is such that, on the periodic orbit, 
the virtual constraints (0 are identically satisfied. How- 
ever, off the periodic orbit, even if the virtual constraints 
are satisfied at the end of given step, they will not in gen- 
eral be satisfied^ at the beginning of the next step. Con- 
sequently, the behavior of the robot cannot be deduced 
from the behavior of the uncontrolled variables q u and the 
simulation of the complete model is required to predict the 
behavior of the robot. In the language of [57], [231 Chap. 5], 
while the feedback control law (fT0|) has created a zero dy- 
namics of the stance phase dynamics, it has not created a 
hybrid zero dynamics, that is, a zero dynamics of the full 
hybrid model (|6|. 

If the control law could be modified so as to create a hy- 
brid zero dynamics, then the study of the swing phase zero 
dynamics (|15p and the impact model would be sufficient to 
determine the stability of the complete closed-loop behav- 
ior of the robot, thereby leading to a reduced-dimension 
stability test. A modification of the control law to achieve 
a hybrid zero dynamics was first proposed in [28] ; a second 
more easily implementable method has been given in |29j . 
along with a complete stability analysis. 

Following [3S] , the virtual constraints are modified stride 
to stride so that they are compatible with the initial state 

7 This may be true for several reasons, one of which is that the 
virtual constraints may not have been chosen to be compatible with 
the impact map. 
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of the robot at the beginning of each step. The new output 
for the feedback control design is 



Vc = h(q, yi, fa) = q a - h d (0) - h c (9, jji, fa). 



(20) 



This output consists of the previous output (0, and a 
correction term h c that depends on j7| evaluated at the 
beginning of the step, specifically, j/j = q a ,i — hd{8i) and 
fa = — ^ggp- Oi, where the subscript "i" denotes the 
initial value for the current step. The values of yi,fa are 
updated at the beginning of each step and held constant 
throughout the step. The function h c is taken to be a 
three-times continuously diffcrcntiablc function of 6 such 
thai© 



h c {9i,yi,fa) 



V, 



(21) 



h c (0,yi,fa) = 0, 



7- 



With h c designed in this way, the initial errors of the output 
and its derivative are smoothly joined to the original virtual 
constraint at the middle of the step. In particular, for 
any initial error, the initial virtual constraint hd is exactly 
satisfied by the end of the step. 

The robot's behavior with the new control law is very 
close to its behavior with the fixed virtual constraint, the 
difference being that the initial output error is zeroed by 
the choice of the function h c instead of being approxi- 
mately zeroed by equation Our choice of h c zeroes 
smoothly the initial error, and the initial virtual constraints 
are joined at the middle of the step. Consequently, the 
torque required to zero the error is less than with a high- 
gain controller, the variation of ground reaction forces is 
reduced, and sliding or take-off is avoided. 

Under the new control law defined by (f20|) . the behavior 
of the robot is completely defined by the impact map and 
the swing phase zero dynamics (|15|) , where hd is replaced by 
hd + h c . The stability of a fixed-point x* can now be tested 
numerically using a restricted Poincare map defined from 
SHZ to SnZ, where Z = {(q,q)\y c (q) = 0, y c (q) = 0} 
and S is the switching surface. The key point is that in 
SnZ, the state of the robot can be represented using only 
three independent variables, x z = [qi,qi,6]'. 

The restricted Poincare map P z : S D Z -> S D Z in- 
duces a discrete-time system x% +1 — P z {x z k ). From [29] . 
for e sufficiently small in (fTO)) . the linearization of P z about 
a fixed-point determines exponential stability of the full- 
order closed-loop robot model. Define b~x z k = xf, — x z *. 
The Poincare map linearized about a fixed-point x z * = 
((%)*! (?<)*> H) gi yes r i se to a linearized system, 



(22) 



where the (3 x 3) square matrix A z is the Jacobian of the 
Poincare map and is computed as follows 



A z = [ A\ A\ A% ] 



3x3 ' 



s In our specific application, we used a fifth order polynomial for 
9i < 9 < ° z ^f ■ continuity of position, velocity and acceleration is 
ensured at = ^'^f , 



where 



A" 



P z (x z * + Ax z ) - P{x z * - Ax z 



2Ax z 



and Axf = 




i = 1,2,3, 

(23) 

The quantities Agi, A^i, 



A9 are small perturbations introduced to calculate the lin- 
earized model, and the denominator in (|23[) must of course 
be interpreted as the scalar perturbation used in computing 
Axf. In general, the calculation of the Jacobian is sensi- 
tive to the amplitude of the perturbation Ax\ . The calcu- 
lation of matrix A z requires six evaluations of the function 
P z (x z * ± Ax z ). Each evaluation of this function is com- 
posed of the reconstruction of the vector x\ from x z with 
x z s SnZ, the calculation of the impact map (HI)-©, the 
calculation of yi,fa, the calculation of h c (0) and the inte- 
gration of the swing phase zero dynamics. A fixed-point 
of the restricted Poincare map is locally exponentially sta- 
ble, if, and only if, the eigenvalues of A z have magnitude 
strictly less than one [23j Chap. 4]. 

B. Example of the periodic motion minimizing integral- 
squared torque 

We consider the virtual constraints corresponding to the 
optimal periodic motion obtained in Section |IV-C[ and the 
control law defined by ([20)) is used. 

To study the stability of this control law around the pe- 
riodic motion, the eigenvalues of the linearized restricted 
Poincare map are computed (for Aqi — 0.0750°, Aqi = 
A8 = O.^s- 1 ), yielding 



A z 



0.1979 
5.8899 
-4.7411 



The three eigenvalues are: 

Ai = 
A 2 = 
A 3 = 



-0.4625 -0.2145 
-2.8417 -1.7476 
-0.2132 0.7809 



0.8878 

-0.6951 

-2.0891 



(24) 



One eigenvalue has magnitude greater than one and hence 
the gait is unstable under this controller. 

We have found that for most periodic motions optimized 
with respect to integral-squared torque per step length, 
(fT9| . the obtained gait is unstable under the control law 
defined by (|2"0"]) . In the following sections, three strate- 
gies to obtain stable walking will be considered. In the 
first strategy, closed-loop stability is directly considered in 
the design of the periodic motion. In the second strategy, 
a stride-to-stride controller is introduced to stabilize the 
gait. In the third strategy, freedom in the selection of the 
controlled outputiH is used to obtain a stable walking cycle 
using only within-stride control. 

9 The controlled outputs are no longer the actuated variables as in 
I10I I. but a judiciously chosen linear combination of q. A convenient 
choice of outputs is given. 
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VI. Periodic motion optimized with respect to a 

STABILITY CRITERION 

Because we are using optimization to compute a periodic 
solution of the model, it is possible to consider the stabil- 
ity condition (i.e., magnitude of eigenvalues less than one) 
either as a criterion in the optimization process, or as a 
constraint. 

Using the stability test above, which has a reasonable 
calculation cost, the selection of a periodic motion that 
minimizes the maximum of the magnitudes of the eigen- 
value of A z was performed, 

J = max{|A(A z )|}, 

where X(A Z ) is the set of eigenvalues of A z . The optimiza- 
tion process is similar to the optimization process described 
in Section ITVl only the criterion is changed. Starting from 
the same initialization for the optimization process, a new 
periodic trajectory is obtained. This trajectory is defined 
by 

q* f = [-0.0306,-0.3304,0.3892,-0.2853,0.0703,0.0265 
-0.4948,0.2827]' 
q* f = [-0.2719, -1.6158, -0.0710, -0.1553, -0.1998, 0.2312 
1.1816,-0.0450]'. 

The stick-diagram for a step is presented in Fig. [101 The 
walking cycle has a period of T = 0.175 seconds, a step 
length of L — 0.144 m, and an average walking speed of 0.82 
m/sec. The step width is 0.173 m, close to the hip width. 
The nominal gait's joint profiles and angular velocities over 
two steps are shown in Figs. [TT] and [12J respectively. The 
unactuated and actuated variables are presented, and it is 
shown that 8 is monotonic. The duration of a step is shorter 
than in the case of torque optimization (i.e., the motion is 
faster). The movement in the frontal plane has smaller 
amplitude; indeed, the maximal sidc-to-sidc variation of 
the center of mass in the frontal direction is less than 1.5 
cm (see Fig. ITU)) , while this deviation was more than 5 cm 
in the case of torque optimization (see Fig. [5]). 
The eigenvalues of the matrix A z are 

Ai, 2 = -0.873 ±216i 
|Ai, 2 | = 0.890 
A 3 = 0.887 

All of the eigenvalues have magnitude less than 1.0, indi- 
cating that the obtained nominal orbit x* is locally expo- 
nentially stable. To illustrate the orbit's local stability, the 
3D biped's full model in closed-loop is simulated with an 
initial state perturbed from the fixed-point x* . An initial 
error of —0.5° is introduced on each joint. 

Fig. [13] shows the evolution of the final values of the 
uncontrolled variables q u from one step to the next. These 
variables converge toward the periodic motion. With the 
modification of the virtual constraints by the introduction 
of the polynomial h Cl the output y c and its derivative are 
zero throughout each step. 




Fig. 10 

Stick-diagram of the optimal motion with respect to the 
stability criterion 




0.1 . 0.2 , 0.4 0.1 . 0.2 , 0.3 0.4 

time (sj time (sj 

Fig. 11 

Joint profiles of the obtained periodic motion over two 
steps. where the small circles represent q*j . 




0.1 . 0.2 . 0.3 0.4 0.1 . 0.2 . 0.3 0.4 

time (sj time (sj 

Fig. 12 

Joint rate profiles of the obtained periodic motion over 
two steps. where the small circles represent . 
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20 30 40 

step number 



T3 
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Cfi 1 .7 a 

-c 

,2,1.69 >■_ 

1.6 1 

1.55 



10 20 30 40 J.„ ._. 

step number step number 

Fig. 13 

The evolution of q u at the end of each step for the 3D 
biped's full model under closed-loop walking control, with 
the initial condition perturbed from x* . the small circles 
represent the values on the periodic orbit. 



Fig. [HI shows the phase-plane evolution of the first four 
variables as an illustration of the behavior of the state of 
the robot. The discontinuity at impact manifests itself as 
a straight line on the plots. The desired motion of the 
robot is periodic, with exchange of the legs at impact. This 
behavior corresponds to a motion with a period of one step 
for the variables q2, (73 and and a period of two steps for 
qi ; see (|TJ) . The convergence towards a periodic motion is 
clear for both the controlled and uncontrolled variables. 





q 2 (rad) 



U.4 -U.H -U.b -U.4 iT -* U 

53 (rad) <J4 (rad) 

Fig. 14 

Phase-plane plots for qt, i = 1, . . . , 4. The straight lines 

CORRESPOND TO THE IMPACT PHASE, WHERE THE STATE OF THE 
ROBOT CHANGES INSTANTANEOUSLY. THE INITIAL STATE IS 
REPRESENTED BY A (RED) STAR. EACH VARIABLE CONVERGES TO A 
PERIODIC MOTION. 



VII. Stride-to-Stride CONTROLLER 

If a desired periodic gait is not stable, or if the corre- 
sponding rate of convergence is not sufficiently rapid, then 



event-based control may be designed and integrated with 
the continuous, stance phase controller [50] . 

A. General method 

Let /3 be a vector of parameters that are held constant 
during the stance phase and updated at each impact. The 
parameters could be a subset of the parameters used to 
specify the virtual constraints, ao, •■■ ,03, or an auxiliary 
set of parameters. Moreover, the parameter updates can 
be done on the basis of the full-state of the robot as in [1] , 
or on the basis of the state of the hybrid zero dynamics, 
which is done here. 

The output in (|20[) is augmented with an additional term 
h s {8,P) depending on a vector of parameters /?: 

V = Hq,Vi,yi,f3) = q a -hd(9)-h c (9,yi,yi)-h s (9,p) (25) 



wit] 



10 



h.(e,p) = 






P 

0, 0.16» 4 



(26) 



0.961/ < 8 < 9 



This choice is convenient because it allows the value of the 
virtual constraints at mid-stance to be updated, without 
requiring a re-design of the event-based controller that cre- 
ated the hybrid zero dynamics. As long as the impact oc- 
curs for 8 > O.Wi + 0.98 f, which happens when the motion 
is close to the periodic orbit, the final state of the robot is 
on the original hybrid zero dynamics independently of any 
modification h s for the previous step. 

The Poincare map can now be viewed as a nonlinear 
control system on S (~l Z with inputs pk, where Pk is the 
value of P during the step k + 1, namely xl +1 — P(x^,Pk)- 
Linearizing this nonlinear system about the fixed point x z * 
and the nominal parameter value P* — 06 xi leads to 



Sx z k+1 = A z Sx z k +F5p k , 



(27) 



where 6Pk = Pk — P* arid F is the Jacobian of P with 
respect to p. Designing a feedback matrix 



Sp k = -K5x z k 



(28) 



such that the eigenvalues of (A z — FK) have magnitude 
strictly less than one will exponentially stabilize the fixed 
point x* . 

B. Example of the periodic motion minimizing the integral- 
squared torque 

The periodic motion described in Section IIV-CI can be 
stabilized using a stride-to-stride controller. The virtual 
constraints are modified according to equation (|25l) ; for 
the fc-th step the virtual constraints hd,k are 



hd,k = h d{8) - h c (8,yi,yi) - h s (9,p k ). 



(29) 



10 In our specific application, we used a sixth order polynomial for 
8i < 6 < O.Wi + 0.90/ . Continuity of position, velocity and accelera- 



tion is ensured at 9 = O.Wi + 0.90 



/■ 
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The (3 x 6) matrix F was computed numerically, analo- 
gously to ([22} , yielding 



-0.030 -0.028 
-0.237 -0.233 
0.163 0.224 



0.141 -0.065 0.028 0.018 
1.073 -0.494 0.225 0.144 
-0.023 0.001 -0.231 -0.064 

(30) 

The matrix A z is given in ([24| . 

The (6 x 3) gain matrix K was calculated via DLQR 
so that the state-feedback law 5/3k = —KSx^. minimized 
the cost function J2k i^ x k ^ x k + r ^P'k^^k), subject to the 
state dynamics (|27|). The coefficient r allows a tradeoff to 
be made between the convergence rate and the amplitude 
of the change of the virtual constraint, which affects the 
region of convergence. For r = 2 we obtain 



K = 



-0.603 


0.243 


0.176 


-0.607 


0.186 


0.171 


2.566 


-1.704 


-0.814 


-1.181 


0.793 


0.376 


0.590 


-0.164 


-0.165 


0.359 


-0.179 


-0.108 



(31) 



The eigenvalues of the linearized Poincare map in closed 
loop A z — FK indicate that the gait is stable with the ad- 
dition of the stride-to-stride controlled 1 *!. The eigenvalues 
become 



Ai 

^2,3 
1-^2,31 



0.7906 

-0.4478 ± 0.047i 
0.4503 



To illustrate the orbit's local stability, the 3D biped's 
complete model in closed-loop is simulated with the initial 
state perturbed away from the fixed-point x* . An initial 
error of —1° is introduced on each joint and a velocity error 
of — 5°s _1 is introduced on each joint velocity. 

Fig. [15] shows the evolution of the values of the uncon- 
trolled variables q u at the end of each step, when event- 
based DLQR controller is used. These variables clearly 
converge toward the periodic motion. 

Fig. [TBI shows phase-plane plots of the first four variables. 
The convergence towards a periodic motion is clear for the 
controlled and uncontrolled variables. 

VIII. Improved Selection of the Outputs to be 
Controlled 

In the previous two sections, the controlled variables 
driven by the virtual constraints are simply the actuated 
variables, q a ; see ([7]). The choice of the controlled vari- 
ables directly affects the zero dynamics in (flor) . It is shown 
here that for the same desired periodic motion, the stabil- 
ity of the closed-loop system can be dramatically improved 
through a judicious choice of the controlled variables. 



n The controller of Sec. IV-BI is used within stride and the event- 
based controller is applied at each impact. 



-o.oi 

-0.02* 
-0.03 
-0.04 



0.19 
0.18' 

"2 0-17 

TO 

' 0.1 6 

'S 
15 15(» 



10 15 20 

step number 



step number 




5 10 15 20 25 5 10 15 20 25 

step number step number 

Fig. 15 

The evolution of g u at the end of each step for the 3D 
biped's full model under closed-loop walking control, with 
the initial condition perturbed from x* . the small circles 
represent the values on the periodic orbit. 





91 (rad) 



<J2 (rad) 




93 T rad ) 




-0.8 -0.6 ,u.- 

<?4 (rad) 



Fie 



16 

Phase-plane plots for qi, i = 1, ... ,4. The straight lines 

CORRESPOND TO THE IMPACT PHASE, WHERE THE STATE OF THE 
ROBOT CHANGES INSTANTANEOUSLY. THE INITIAL STATE IS 
REPRESENTED BY A (RED) STAR. EACH VARIABLE CONVERGES TO A 
PERIODIC MOTION. 



A. Effect on the swing phase zero dynamics 

For simplicity, we limit our analysis to the case of con- 
trolled variables that are linear with respect to the config- 
uration variables. Thus the controlled variables are 



M 





q„ 



[ M 1 M e M a } 



0.1 



q a 



(32) 



where M is a (6 x 8) constant matrix with M a invertible. 
A known periodic motion q*(t) can be reparameterizecQ 
as function of the variable 0, yielding q*(0). The virtual 



12 This assumes that 9 is monotonic. 
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constraint for the new controlled variables then yields the 
output 





9i 








y = M 


6 


-M 


9 


(33) 












ho(q) 




h d {6) 





When the constraint is satisfied, y = 0, equation 
allows us to solve for q a , giving 

q a = q* a (6) + M- 1 M 1 (ql(6)-q 1 ). (34) 

Substituting this equation into (I14[) , we obtain for the 
swing phase zero dynamics 



Du(q u 



D 12 {q u )( 9 -^t 



Ojh. 

d9 



+ M 2 - q \)+H 1 (q u ,q. 



= 0. 

(35) 

The nominal periodic motion satisfies both equations 
(|35|) and (p~5|) . but the two equations produce different so- 
lutions away from the periodic motion. When the principle 
of virtual constraints is applied to a system with only one 
degree of underactuation, namely 9, which is common for 
example in planar bipeds, the swing phase zero dynamic is 
not affected by the choice of the output, and therefore the 
stability of a periodic orbit (i.e., walking motion) is not 
modified; only the transient motion can be different. In 
other words, for planar robots with one degree of underac- 
tuation, the stability depends only on the trajectory of the 
periodic orbit and not on the choice of virtual constraints 
used to achieve it [31], [23l pp. 160]. 

In the case of a system with two degrees of underac- 
tuation, the choice of the controlled output can affect the 
stability of the gait via the choice of M~ X M\ . In order to il- 
lustrate this property, a new choice of output is proposed. 
This choice is based on the following physical reasoning: 
The motion in the frontal direction is difficult to stabilize. 
The position of the center of mass in the frontal direction 
is important. If, at touchdown, the center of mass is not 
between the feet, but outside the position of the next sup- 
porting foot, the robot will topple sideways. Thus, the 
control of the variable q^ (which regulates step width on 
the swing leg) is replaced by the control of the distance be- 
tween the swing leg end and the center of mass along the 
frontal direction. To obtain a linear output, this function 
is linearized around the touchdown configuration to define 
M in (1331). 



B. Example of the periodic motion minimizing integral- 
squared torque 

The periodic motion described in Section IIV-CI can be 
stabilized using the new controlled output. As mentioned 
in the previous subsection, the actuated joints (73, (74, (75, (77, 
and q% are controlled via virtual constraints just as in the 
original control law. A new output hdA is considered, with 
this output no longer based on q e but instead on distance 
between the swing leg end and the center of mass along the 
frontal direction. 



For this trajectory, for support on leg 1, the linearization 
around qf of the distance between the swing leg end and 
the center of mass along the frontal direction yields 

d = -0.457qi - 0.020q2 - 0.018g 3 - 0.020*7 4 - 0.489g 5 
+0.461*76 - 0.056g 7 - 0.022g 8 . 

On the periodic orbit, this distance is evaluated and ap- 
proximated by a function of 9, denoted d*{9). The new 
controlled output is then 

y 4 = -0.457gi - 0.020<7 2 - 0.018*7 3 - 0.020*? 4 - 0.489q 5 
+0.461*7 6 - 0.056*77 - 0.022g 8 - d*(9). 

(36) 

When the control law is defined using this new output, 
the walking gait is stable, as can be shown via the calcula- 
tion of the eigenvalues of the linearization of the restricted 
Poincare map, A z . The eigenvalues are: 



Ai 

^2,3 
1-^2,3 1 



0.7846 

-0.028 ±0.250i 
0.2512. 



To illustrate the orbit's local exponential stability, the 3D 
biped's model in closed-loop is simulated with the initial 
state perturbed from the fixed-point x*. An initial error 
of —1° is introduced on each joint and a velocity error of 
— 5°s _1 is introduced on each joint velocity. 

Fig. [T7] shows the evolution of values of the uncontrolled 
variables q u at the end of each step when the new output is 
used. These variables clearly converge toward the periodic 
motion. 
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Fig. 17 

The evolution of q u at the end of each step for the 3D 
biped's full model under closed-loop walking control, with 
the initial condition perturbed from x* . the small circles 
represent the values on the periodic orbit. 



Fig. [181 shows phase-plane plots of the first four variables. 
The convergence towards a periodic motion is clear for the 
controlled and uncontrolled variables. 
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Fig. 18 

Phase-plane plots for. g i} i = 1, . . . ,4. The straight lines 

CORRESPOND TO THE IMPACT PHASE, WHERE THE STATE OF THE 
ROBOT CHANGES INSTANTANEOUSLY. THE INITIAL STATE IS 
REPRESENTED BY A (RED) STAR. EACH VARIABLE CONVERGES TO A 
PERIODIC MOTION. 



IX. Conclusions 

A simple 3D bipedal model has been studied, with the 
objective of developing a time-invariant feedback control 
law that induces asymptotically stable walking, without 
relying on the use of large feet. For this reason, a biped 
consisting of five links, connected to form two legs with 
knees and a torso, was assumed to have point feet with 
no actuation between the feet and ground. Inspired by 
its success in solving similar problems for planar robots, 
the method of virtual constraints was applied to the 3D 
robot, with the virtual constraints chosen via optimization 
as suggested in [5]. 

The main contributions of the paper are: 

1. The development of an efficient optimization process to 
obtain periodic motions in 3D for a robot with point con- 
tact feet. 

2. The computation of human- like periodic walking mo- 
tions that can be stable or unstable, depending on the 
choice of actuated variables and corresponding virtual con- 
straints. 

3. The numerical study of stability on the basis of a low- 
dimensional subsystem corresponding to the hybrid zero 
dynamics. The Poincare return map was computed in a 
space of dimension three for a robot with two degrees of 
underactuation. 

4. The use of a stride-to-stride controller to stabilize a 
walking motion that was not naturally stable for a given 
choice of controlled outputs. 

5. The discovery of the importance of the selection of the 
controlled outputs on the stability of a given periodic mo- 
tion. 

Points (1), (3) and (4) can be viewed as direct exten- 
sions of the work previously done on planar walking robots. 
However, the points (2) and (5) are specific to the study of 
robots in 3D. 



The control strategy developed here for the control of 
robots with point feet can be extended to the case of robots 
with actuated, non-trivial feet, as was done in |17j . |32j . 
[2"2"] for planar robots. In this case, foot rotation can be 
included as part of a natural gait and the explicit control 
of the ZMP position can be addressed. Another possible 
extension would be to continue with the unactuated point- 
foot hypothesis, this time including yaw rotation about the 
foot. It is unclear at this time whether this extension is 
straightforward or not. 
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